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Computer  program  RMA-7  (1(+n^etTrm973f  has  been  expanded  to  be  able 
nl!!ulate  den?1ty  *"duced  and  water  quality  conditions  typically  found 

Iwin^i?  rS *.i  T*e  mathemat1cal  bas1$  for  the  program  and  methods  of 
implementing  the  code  for  various  kinds  of  flow  are  discussed.  Results  from 

di.ffere"t  applications  are  presented.  The  first  example  simulates  flow 
conditions  that  were  measured  in  a  physical  model  of  a  deep  reservoir.  Results 
are  also  Included  for  simulations  of  the  circulation,  temperature  (CONTINUED) 


,20  (CONTINUED) 

and  dissolved  oxygen  concentrations  for  Lake  Taneycomo,  Missouri:  Comparison 
of  measured  and  computed  results  from  both  applications  shows  that  the 
RMA-7  model  provides  reasonably  accurate  results  for  both  flow  circulating 
and  water  quality. 
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APPLICATION  OF  THE  FINITE  ELEMENT  METHOD  TO  VERTICALLY 
STRATIFIED  HYDRODYNAMIC  FLOW  AND  WATER  QUALITY1 

2  3 

R.  C.  MacArthur  and  W.  R.  Norton 


INTRODUCTION 

Continuing  interest  in  the  internal  processes  of  reservoirs,  lakes  and 
estuaries  has  intensified  the  development  of  mathematical  models  for 
simulation  of  vertically  stratified  flow.  Motivated  primarily  by  the 
long  term  desire  to  predict  the  water  quality  response  to  system 
modifications,  current  modeling  efforts  have  focused  on  the  need  to 
describe  flow  coupled  with  temperature  and/or  salinity  fields  in  order 
to  forecast  the  influence  of  density  induced  flows*  Computational 
algorithms  have  shown  sufficient  promise  that  efforts  are  under  way 
to  collect  prototype  data  which  can  be  used  for  calibration  and 
verification  of  both  flow  and  water  quality  models. 

Among  the  several  models  which  have  been  formulated  to  simulate  density 
induced  flows  and  water  quality  is  one  called  RMA-7.  This  model,  which 
was  originally  developed  for  the  Office  of  Water  Resources  Research, 

(King,  1973)  has  existed  for  several  years  but  has  received  relatively 
little  use  in  prototype  applications.  It  is  the  intent  of  the  work 
reported  herein  to  demonstrate  the  operation  of  the  current  version  of 
RMA-7  and  to  offer  appropriate  information  and  comments  on  the  use  and 
implementation  of  the  model. 

Specifically,  the  paragraphs  which  follow  contain  a  brief  description  of 
the  mathematical  basis  of  RMA-7  as  well  as  an  example  of  its  implementation 
on  a  prototype  system  which  approximates  the  physical  dimensions  of  the 
so-called  GRH  flume  at  the  Corps’  Waterways  Experiment  Station  in  Vicksburg. 
Also  included  are  some  results  obtained  by  the  model  from  the  simulation  of 
temperature  and  dissolved  oxygen  for  Lake  Taneycomo  in  Missouri.  Several 
statistical  comparisons  are  included  between  simulated  and  measured  values 
for  Lake  Taneycomo  which  are  designed  to  quantify,  to  some  extent,  the 
accuracy  of  the  model. 


^Presented  at  the  3rd  International  Conference  on  Finite  Elements  in  Water 
Resources,  May  19-23,  1980,  the  University  of  Mississippi,  Oxford,  MS. 
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GOVERNING  EQUATIONS  -  DIFFERENTIAL  FORMS 


RMA-7  is  a  two  dimensional  mathematical  model  which  describes  the  behavior 
of  velocity,  pressure,  temperature  and  dissolved  oxygen  in  the  vertical 
plane  with  homogeneity  assumed  in  the  third  (lateral)  direction;  the  model 
will  accommodate  width  gradients  in  both  the  X  and  Y  directions.  This 
model  is  formulated  on  the  classical  concepts  of  conservation  of  mass, 
momentum,  and  energy  although  it  is  somewhat  unusual  as  it  retains  the 
complete  vertical  momentum  equation  and  does  not  make  the  assumption  that 
pressure  must  be  hydrostatic. 


Hydrodynamic  Model 

The  equations  used  in  RMA-7  have  the  following  differential 
forms: 


Velocity  equations. 
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Continuity  equation. 
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Temperature  equation. 
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The  equations  used  in  RMA-7  to  describe  the  behavior  of  dissolved  oxygen 
have  the  differential  forms: 


Dissolved  oxygen  equation. 
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demand  (BOD)  equation 

Phytoplankton  equation 


where  C\  •  dissolved  oxygen  concentration 
c 2  *  carbonaceous  BOD  concentration 
Cj  •  phytoplankton  concentration  (dry  weight) 

K2  *  BOD  decay  rate 

a  -  the  mass  of  oxygen  produced  per  unit  mass  of  phyto¬ 
plankton  growth 

V  “the  local  rate  of  phytoplankton  growth 
“  Vmax  (I  +  Ki*  "  r 

where  Umay  "  the  maximum  specific  phytoplankton  growth  rate  at 
the  local  temperature 

I  “  the  intensity  of  light  at  the  local  depth 
Kx  “  the  light  half  saturation  constant 
r  “  the  phytoplankton  respiration  rate  at  the  local 
temperature 

and  all  other  terms  as  previously  defined. 

Equations  1  through  7  represent  the  mathematical  basis  for  the  model  RMA-7. 
It  is  recognized  that  these  equations  provide  an  approximate  representation 
of  the  governing  processes  in  that  certain  second  order  terms  with  respect 
to  width  have  been  dropped  and  that  the  viscosity/dispersion  terms  are 
largely  empirical.  Notwithstanding  these  shortcomings,  however,  the  above 
relationships  have  shown  promise  in  simulation  of  observed  phenomena,  and 
provide  a  general  framework  for  testing  and  improving  the  procedures 
necessary  for  simulation  of  vertically  stratified  flow. 

As  can  be  seen,  the  first  four  equations  and  the  second  three  equations 
each  form  a  closed  set  which  must  be  solved  simultaneously.  The  coupling 
between  the  two  sets  is  manifest  in  the  convective  velocities  u  and  v, 
and  in  the  water  temperature,  T.  Fortunately,  the  coupling  between  the 


two  sets  is  of  the  nfeed  forward"  type  in  that  equations  1  through  4  may 
be  solved  independently  of  equations  5,  6  and  7,  with  the  results  of  the 
hydrodynamic  solution  acting  like  coefficients  in  the  solution  of  the 
water  quality  model. 

In  addition  to  the  volumetric  terms  shown  above,  the  model  contains  a 

number  of  additional  features  to  account  for  surface  effects  and  source/ 

sink  terms.  Host  important  of  these  are:  1)  bottom  friction  at  the  f 

soil-water  interface  calculated  as  a  function  of  bottom  velocity  and  a 

Chezy  coefficient;  2)  surface  heat  exchange  as  a  function  of  the  water 

temperature,  the  equilibrium  temperature  and  an  exchange  coefficient; 

3)  internal  absorption  of  short  wave  solar  radiation  as  a  function  of 
depth;  and  4)  surface  oxygen  exchange  as  a  function  of  the  local  oxygen 
deficit  and  an  exchange  coefficient  which  is  a  function  of  wind  speed. 

SOLUTION  TECHNIQUE 

The  governing  equations  are  solved  by  the  finite  element  method  using 
Galerkin's  criteria  for  the  method  of  weighted  residuals.  The  formulation 
employs  a  mixed  set  of  basis  functions,  with  quadratic  functions  used  for 
all  state  variables  except  pressure  where  a  linear  function  is  used.  The 
linear  pressure  function  implies  a  constant  element  density,  which  is 
calculated  as  a  function  of  average  nodal  temperatures.  Greenfs  Theorem 
is  used  to  lower  the  order  of  all  second  derivatives  in  the  viscosity/ 
dispersion  terms,  resulting  in  surface  integrals  which  must  be  evaluated 
(either  implicitly  or  explicitly)  along  system  extremities.  Green's 
Theorem  is  also  used  on  the  pressure  terms  of  equations  1  and  2  permitting 
a  surface  integral  to  be  used  as  a  discharge  boundary  condition  rather  than 
a  specified  pressure  value;  this  procedure  permits  retention  of  all  nodal 
continuity  equations  and  substantially  improves  the  model's  performance. 

RMA-7  uses  an  implicit,  Newton-Raphson  computation  scheme  to  achieve  a 
solution  to  the  set  of  nonlinear  equations  which  define  the  model.  The 
computer  program  accommodates  either  triangular  or  quadrilateral 
isoparametric  elements  with  numerical  integration  used  to  evaluate 
all  surface  and  area  integrals.  The  use  of  the  isoparametric  formulation 
with  interelement  geometric  slope  continuity  allows  flow  to  move  parallel 
to  boundaries  at  all  points  while  at  the  same  time  providing  a  means  for 
reasonable  representation  of  real  physical  systems. 

EXAMPLE  PROBLEMS 

The  model  described  above  has  been  applied  to  several  physical  situations, 
two  of  which  are  summarized  below.  The  first  example  shows  results 

calculated  from  tests  conducted  with  the  geometry  of  the  GRH  flume  at  WES.  * 

The  second  shows  results  calculated  from  conditions  found  in  Lake  Taneycomo 
in  Missouri.  The  first  example  was  run  using  only  the  flow  and  temperature 
portions  of  the  model,  while  the  second  includes  flow,  temperature  and 
dissolved  oxygen.  The  networks  used  for  each  of  these  example  problems 
are  reproduced  in  figure  1. 
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Example  Problem  1  -  GRH  Flume 


The  GRH  flume  is  a  hydrodynamic  test  facility  which  is  80  feet  long  and 
varies  in  depth  from  approximately  1  foot  at  its  upper  (inflow)  end  to 
3  feet  at  its  lower  (outflow)  end,  with  two  different  cross  sections 
along  its  length*  The  inflow  section,  which  is  20  feet  long,  has  a 
horizontal  bottom  and  varies  linearly  in  width  from  1.0  foot  to  2.85 
feet.  The  lower  section  has  a  constant  width  of  2.85  feet,  but  a  bottom 
which  drops  2  feet  over  its  60  foot  length. 

RMA-7  was  applied  to  the  GRH  flume  geometry  by  construction  of  a  finite 
element  network  containing  57  elements  and  158  node  points  as  shown  in 
figure  1.  In  constructing  this  network  it  was  felt  desirable  to  allow 
flow  to  move  parallel  to  the  flume  bottom  at  all  locations.  For  this 
reason  continuous  curves  were  passed  through  the  breakpoint  on  the  flume* s 
bottom  20  feet  from  the  upstream  end  and  at  the  transition  to  the  outlet. 
This  type  of  construction  permits  continuous  velocities  to  exist  along 
the  flume  bottom  and  completely  eliminates  artificial  stagnation  points. 
The  small  discharge  nozzle  at  the  outlet  has  been  included  for  each 
boundary  condition  specification  and  does  not  exist  on  the  physical  flume. 

To  run  RMA-7  it  is  necessary  to  specify  values  for  eddy  viscosity  and 
eddy  diffusion  coefficients.  At  present,  this  is  done  by  experience 
and  numerical  testing  of  problems  which  have  known  or  assumed  velocity 
and  temperature  distributions.  In  the  GRH  flume  case,  a  series  of 
numerical  tests  were  conducted  on  a  steady  state  problem  to  determine 
a  set  of  satisfactory  coefficients,  and  the  network’s  sensitivity  to 
the  various  coefficients.  The  values  determined  for  the  GRH  flume  had 
the  relative  values  of  exx  =  0.05,  eXy  =  0.0005,  eyX  =  0.01,  Cyy  =  0.1, 

Dx  *  0.1,  and  Dy  =  0.0005,  after  accounting  for  element  distortion  and 
size. 

Results  from  two  examples  are  shown,  one  for  a  homogeneous  flow  and  one 
for  a  nonhomogeneous  flow.  In  each  case  a  flow  of  10  gpm  was  introduced 
into  a  still  flume  with  a  linear  velocity  distribution  in  the  lower 
element  at  the  inflow  end.  The  homogeneous  case  was  run  isothermally 
at  a  temperature  of  10.3QC,  while  the  nonhomogeneous  case  was  started 
with  an  initial  temperature  of  10.3°C  and  an  inflow  of  5°C. 

Velocity  distributions  produced  by  each  condition  are  graphically  compared 
at  three  times  in  figure  2,  The  effects  of  the  density  stratification  are 
quite  evident  between  the  two  cases  with  the  colder,  more  dense  water 
underflowing  the  lighter  and  warmer  water  near  the  surface.  The  results 
shown  are  representative  of  ongoing  work  with  the  GRH  flume,  although  no 
measured  data  is  currently  available  for  model/prototype  comparisons 
under  the  conditions  simulated.  The  general  shape  of  velocity  distribu¬ 
tions,  and  the  arrival  time  of  the  temperature  front  (15-18  min)  are  in 
general  agreement  with  measured  data,  however,  and  suggest  the  model  will 
perform  well  when  suitable  data  become  available. 
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STATION  2 _ STATION  1 


Figure  1. — Example  Problem  Networks 


Figure  2. — Comparison  of  Homogeneous  and  Non  Homogeneous  Simulated  Velocity 
Distributions  for  the  GRH  Flume  Example  Problem 


K xamp 1 e  V rob  Lem  2  -  Liko  Taney  coim > 

The  second  example  problem  presents  results  obtained  from  the  simulation 
of  Lake  Taneycomo  in  southern  Missouri.  Lake  Taneycomo  is  a  26-mile-long 
reservoir  bounded  by  Table  Rock  Dam  upstream  and  Bull  Shoals  Reservoir 
downstream,  and  is  used  for  power  production  and  recreation  among  other 
things.  As  low  dissolved  oxygen  has  been  observed  in  the  lake,  RMA-7  was 
applied  with  the  goal  of  evaluating  the  impact  of  changes  in  reservoir 
operation  on  the  ambient  levels  of  dissolved  oxygen. 

To  conduct  the  required  simulations  a  network  of  92  elements  and  246  node 
points  was  constructed  as  shown  in  figure  1;  lateral  width  varied  from 
about  200  to  1000  feet  with  depth  at  a  typical  cross  section.  Diffusion 
coefficients  and  eddy  viscosity  coefficients  were  chosen  to  be  consistent 
with  those  used  in  the  GRH  flume  when  scaled  for  element  distortion. 

Detailed  water  temperature  and  dissolved  oxygen  measurements  were  available 
for  three  separate  week  long  periods  in  the  fall  of  1977.  RMA-7  was 
calibrated  against  two  of  these  periods  and  verified  against  the  third. 
Typical  measured  and  simulated  vertical  profiles  for  temperature  and 
dissolved  oxygen  are  presented  for  stations  T6  and  T8  in  figure  3  for 
the  earliest  calibration  period.  It  should  be  noted  that  the  model  had 
simulated  over  five  days  of  operation  by  the  time  shown  in  these  figures, 
and  that  the  inflow  hydrograph  varied  from  0  to  7000  cfs  in  a  four  to  six 
hour  period  on  a  regular  basis. 

In  order  to  bring  a  degree  of  quantification  to  the  accuracy  of  the  model, 
a  linear  least  squares  regression  of  simulated  to  observed  temperature 
and  dissolved  oxygen  has  been  made  as  shown  in  figure  4,  with  the  statistics 
for  each  fit  given  in  table  1. 

As  can  be  seen,  these  statistics  indicate  a  fairly  good  fit  of  both  temper¬ 
ature  and  dissolved  oxygen,  and  seem  reasonable  in  light  of  the  uncertainty 
in  both  the  model* s  upstream  inputs  (BOD,  dissolved  oxygen,  meteorological 
data,  etc.)  and  the  usual  measurement  errors. 


Table  1. — Regression  Statistics 


Regression 

Statistic 

STATION 

T6 

STATION 

T8 

Temperature 

Oxygen 

Temperature 

Oxygei 

intercept 

4.77 

3.27 

0.29 

3.30 

slope 

0.58 

0.43 

0.96 

0.39 

s.d.  error 

0.65 

0.49 

0.43 

0.30 

correlation 

0.44 

0.30 

0.89 

0.81 
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TEMP  DEG  C  *10’ 

Station  T6 


TEMP  DEG  C  *10 

Station  T8 


(MG/U  SIM  TEMP  (DEG  C) 


SUMMARY  AND  CONCLUSIONS 


The  above  information  summarizes  the  application  of  a  two  dimensional 
finite  element  model  (RMA-7)  to  two  prototype  stratified  flow  situations. 
In  each  application  a  stratified  flow  was  simulated  as  a  result  of  density 
differences  arising  from  temperature  gradients.  In  the  second  example 
dissolved  oxygen,  BOD  and  phytoplankton  were  also  routed  in  accordance 
with  the  overall  flow  fields. 

Based  on  the  data  contained  herein  it  seems  reasonable  to  conclude: 

.  the  finite  element  method  in  general  and  RMA-7  in  particular, 
is  capable  of  simulating  vertically  stratified  two  dimensional 
flow. 

.  the  proper  definition  of  eddy  viscosity  and  dispersion  coeffi¬ 
cients  is  essential  to  proper  model  operation,  and  that  there 
may  be  transferability  of  values  from  problem  to  problem  if 
element  size  and  distortion  is  accounted  for. 

.  the  RMA-7  model  appears  to  give  reasonable  answers,  but  it  is 
not  currently  possible  to  evaluate  its  accuracy  in  relation  to 
velocity  due  to  a  lack  of  prototype  data;  initial  comparisons 
of  temperature  and  dissolved  oxygen  data  are  promising  and 
will  improve  as  calibrations  become  more  rigorous  and 
quantified. 
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